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Abstract. Turbulent motions in the interior of a star play an important role in its evolution, since 
they transport chemical species, thermal energy and angular momentum. Our overall goal is to 
, construct a practical turbulent closure model for convective transport that can be used in a multi- 

dimensional stellar evolution calculation including the effects of rotation, shear and magnetic fields. 
Here, we focus on the first step of this task: capturing the well-known transition from radiative heat 
transport to turbulent convection with and without rotation, as well as the asymptotic relationship 
between turbulent and radiative transport in the limit of large Rayleigh number. We extend the 
closure model developed by Ogilvie (2003) and Garaud and Ogilvie (2005) to include heat transport 
I i and compare it with experimental results of Rayleigh-Benard convection. 
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INTRODUCTION 



> 

■ Turbulent convection plays an essential role in the evolution of most stars. Turbulent 

! motions are often the dominant mechanism of energy transport in stellar convective 

O I zones, and can increase momentum transport and chemical mixing by several orders 

0> ■ of magnitude. In many rotating astronomical systems where turbulence is anisotropic, 

^ . Reynolds stresses are the dominant transporters of angular momentum and therefore 

I influence the internal dynamics of the whole system. 

It is currently not possible to perform a 3-D numerical simulation of convective mo- 
tions over the evolution timescale of a star. However, when considering stellar evolution 
we are not necessarily interested in the specific details of the convective motions, but 
[ rather in statistical properties such as the convective flux and Reynolds stresses. Our 

long-term goal is to construct a closed set of evolution equations for these statistical 
quantities in terms of large-scale system properties (e.g. viscosity, rotation, shear), which 
can be used in a stellar evolution calculation. 

In this paper we focus on modelling how the convective Reynolds stresses and heat 
fluxes are affected by rotation. We are especially interested in approximately predicting 
the onset of convection as well as the asymptotic behaviour of turbulent transport for 
large Rayleigh number; capturing the onset is required for the stellar evolution model 
to correctly place the boundary between the convective and radiative zones. Adequate 
prediction of the asymptotic behaviour of the convective turbulence is important to 
describe the amount of energy and angular momentum that is transported in the majority 
of the convective zone. 

To quantify the quality of the proposed closure model we test it against linear stability 



A practical model of convective dynamics for stellar evolution calculations 



February 1,2008 



1 



analysis, numerical simulations and laboratory experiments of the rotating Rayleigh- 
Benard problem. 

ROTATING RAYLEIGH BENARD CONVECTION 

The typical Rayleigh-Benard convection setup is as follows: two rigid, ideally infinite, 
horizontal plates separated by a distance D confine a weakly compressible fluid between 
them. A temperature difference AT is maintained between a hot bottom and a cool top. 
The system is assumed to be rotating with average angular velocity Q. = Q.z with gravity 
g = —gz. The Boussinesq approximation is valid in this system. The governing equations 
are 

dm = 0, (1) 

{dt + Ukdk)ui + 2£ijkajUk = -a@gi-di^ + vdkkUi, (2) 
idt + Ukdk)@ = Kdkk® (3) 

where the dynamical variables are the temperature offset 0, the pressure perturbation 
from hydrostatic equilibrium ^ and the flow velocity u. The following parameters are 
assumed to be constant: the coefficient of expansion a, the kinetic viscosity V and the 
thermal diffusivity K. Sums over repeated indices are implied. 

The qualitative behaviour of the system is controlled by three dimensionless quanti- 
ties: the Rayleigh number, Ra = agATD^ / (kv) (measuring the ratio between buoyancy 
forces to viscous stabilising forces), the Taylor number, Ta = 4D^Q?/ (measuring the 
ratio between centrifugal forces to viscous forces), and the Prandtl number, Pr = v/k 
(measuring the ratio between the viscous diffusion rate and the thermal diffusion rate). 
For any given Taylor number and Prandtl number, there exists a critical Rayleigh number 
(Rac) above which the system is convective and below which the system is conductive. 
We strive to construct our model to match the known variation of the critical Rayleigh 
number with Taylor number and Prandtl number. We also verify that at Rayleigh num- 
bers much larger than critical the Nusselt number Nu (ratio of total heat flux to conduc- 
tive heat flux) is correctly predicted by the model. 

CLOSURE MODEL 

Our goal is to capture the behaviour of turbulent convection in a rotating system both in 
terms of the onset of turbulence and of its asymptotic properties. We developed a second 
order closure model using the technique described by Ogilvie (2003) and Garaud and 
Ogilvie (2005). 

We write each quantity as the sum of a mean and fluctuating part (u = u + u', ^ = 
+ M^', and = + 0') - in the Rayleigh-Benard problem, these mean quantities only 
vary with z. 

We define the following correlation quantities: Rij = u[u'j, Fi = Q'u'^, and Q = &&, 
so that R = Rii is twice the mean turbulent kinetic energy. 
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The exact equations governing the mean quantities are 

diUi = (4) 

(dt + Ukdk) Ui + lEijk^jUk = - a@gi - diV + vdkkUi - djRij (5) 

{dt + Ukdk)@ = Kdii® - dkFk (6) 

while those governing the second order correlation terms are modelled as 

{dt + Ukdk)Rij + Rikdku j + Rjkdm + lEiiin^iR jm + 2e jimQiRim + a {Fig j + Fjgi) 

-vdkkRij^-CiT-'Rij-C2Z-\Rij-^R5ij)-vCyL-%j, (7) 

{dt + Ujdj)Fi + Rijdj® + FjdjUi + lEijk^jFk +aQgi-h^V+K) dkkFi 

= -CeX-^Fi -\{v + K)CyFiL-^, and (8) 
{dt + Uidi)Q + 2Fidi& - KdkkQ = -CjT-^Q- KCyQL-^ (9) 

where the left-hand- side of each equation is exact, while the right-hand-side models the 
effect of higher order correlations. The constants Ci,C2,Cv,C(,,C-] are free parameters 
of the closure model^ The variable T is the characteristic timescale for the redistribution 
of energy along the turbulent cascade, which is controlled by the turnover time of the 
largest eddies d^^R^/^ and L is the characteristic size of the perturbations near onset, 
which we now describe in more detail. 

Following Prandtl's mixing length theory, Garaud and Ogilvie (2005) originally sug- 
gested that in a wall-bounded experiment L can be thought of as the distance to the wall. 
However, rotation does no net work on the system when Q. is parallel to g and since the 
closure model is constructed on energetic arguments, using this lengthscale prescription 
here fails to capture the known effects of rotation on the onset and turbulent properties 
of convection (Chandrasekhar (1961)). 

Physically, rotation influences the onset of convection by decreasing the characteristic 
lengthscale of convective motions in the direction perpendicular to the rotation axis. To 
capture this effect, we construct our lengthscale L to be the harmonic mean between the 
distance to the wall and the wavelength of the most linearly unstable mode X which, in 
the Rayleigh-Benard experiment, is a function of the Taylor number: therefore 

/I 2\-^'^ 
RESULTS 

We seek solutions to equations (4) - (10) assuming no-slip boundary conditions and 
fixed uniform plate temperatures. The no-slip boundary condition states that ui = at 



Note that Ci ~ 0.4, C2 — 0.6, and Cy — 12 have already been found to give an adequate description of 
the turbulent stresses in Couette-Taylor experiments by Garaud and Ogilvie (2005), but Cg and C7 remain 
to be determined. 
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FIGURE 1. The Rayleigh-Nusselt relation adapted from Chandrasekhar (1961). The dashed line is the 
model prediction for i2 = 0, V = 10~^, Jf = 10~^ generated by varying the temperature difference between 
the two plates. 



the boundaries, implying Rij — Fi = 0. At the lower boundary = AT and at the upper 
boundary = 0. Since the temperature perturbations are zero at both boundaries, Q = 0. 

In Figure 1 we show the Nusselt number - Rayleigh number relationship in the non- 
rotating case for both our model and the selection of experimental data summarised by 
Chandrasekhar (1961). Note that our model shows good agreement with the experimen- 
tal data for the critical Rayleigh number where the transition between conductive (Nu = 
1) and convective (Nu > 1) heat transport occurs. Our model also reproduces the stan- 
dard power law relationship Nu oc RaV3 at high Rayleigh number which is a natural 
consequence of our selection Loc d when Cl = (cf. Prandtl's mixing length theory). 
However, the good quantitative agreement between the experimental data and the model 
prediction for Nu was unexpected since Cg or C7, which play a role when the system 
is convective, have not yet been adjusted from their default value of unity. By adjusting 
these we should be able to further improve the correspondence of the model to reality 
for Ra > Rac. 

In rotating systems, convection is known to delay of the onset of convection at 
high Taylor numbers. Linear theory (cf. Chandrasekhar (1961)) predicts that Rac °« 
Ta^/^, with a coefficient of proportionality which depends somewhat on the boundary 
conditions and on the Prandtl number of the system. In Figure 2, we compare this known 
relationship to our model (solid line). The predicted power law matches linear theory, 
but the coefficient of proportionality is somewhat smaller than required. Nonetheless, 
we consider the agreement satisfactory considering the simplicity of this closure model. 

To conclude, we find that the closure model adequately describes both the turbulent 
convective heat flux as a function of Rayleigh number in the absence of rotation as well 
as the delay of the onset of convection in a rotating system. In future work, we intend to 
compare the degree to which the model matches numerical simulations of developed 
convection in the rotating Rayleigh-Benard problem (e.g. Julien et al. (1996)). This 
comparison will be useful in determining the utility of our model away from the onset 
of convection in a rotating system and permit the calibration of and C7, the remaining 
free parameters of the model. 
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FIGURE 2. The Taylor-Rayleigh critical relation. The solid line is the model prediction for the critical 



Rayleigh number as a function of Taylor number with v = 10 ^, JC = 10 The dotted lines are the 
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asymptotic relations derived from linear stability analysis by Chandrasekhar (1961) for i2 : 
Ta — > oo for the direct mode of instabiUty (see equation 184 on page 106). 



and for 
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